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The thermodynamic and kinetics of the one dimensional lattice gas with repulsive interaction is 
investigated using transfer matrix technique and Monte Carlo simulations. This simple model is 
shown to exhibit waterlike anomalies in density, thermal expansion coefficient and self diffusion. An 
unified description for the thermodynamic anomalies in this model is achieved based on the ground 
state residual entropy which appears in the model due to mixing entropy in a ground state phase 
transition. 
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I. INTRODUCTION 

Water is special fluid for its biological relevance and 
technological applications but most intriguing is that it 
presents thermodynamic and dynamic properties proper¬ 
ties with anomalous (or unusual) behavior [ij. The origin 
of its anomalous properties is actively discussed in the 
literature, with different thermodynamic scenarios com¬ 
peting to describe its behavior on regular and metastable 
regimes Hi- Among alternative views on water ther¬ 
modynamics it should be relevant to mention the second 
critical point hypothesis Q and the singularity free sce¬ 
nario (5|, @ , which will be relevant in the context of the 
current work. 

The second critical point hypothesis was first proposed 
by Poole et al. based on computer simulations of atom- 
istically detailed water models, and relates the observed 
divergence on water’s thermodynamic response functions 
to the critical end point of a deeply metastable liquid- 
liquid phase transition Q|. Singularity free scenario was 
proposed by Sastry et al. using lattice models for liquid 
water which supported the idea that waterlike anomalies 
indeed exist, at least in lattice systems, without a liquid- 
liquid coexistence [§J. In some simplified models of fluid 
it was possible to obtain both second critical point and 
singularity free scenario by adjusting some physical pa¬ 
rameters of the system [3, 0 For instance, in the model 
proposed by Franzese et al. Q it is possible to observe the 
liquid-liquid coexistence to shrink and the temperature 
of second critical point decrease to T c — > 0 by decreas¬ 
ing the correlation between hydrogen bonds inside within 
individual molecules. 

Lattice models of fluid have been extensively used to 
investigate the above mentioned anomalous properties of 
water due to possibility of obtaining analytical or numer- 
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ical results, while exploring a wide range of physical pa¬ 
rameters. In this direction, both thermodynamics |5[ 10] 
and kinetics El III where investigated in lattice models 
with waterlike behavior. Nevertheless, approximations 
employed in two and three dimensions(3D), and even 
some exact solutions in one dimension (ID), tend to gen¬ 
erate complex sets of equations whose analyses is often 
performed numerically. Thus, it should desired to design 
models for which one could obtain simple analytical ex¬ 
pressions connecting thermodynamic anomalous behav¬ 
ior to phase transitions and critical behavior. 


To achieved this goal we previously investigate ID lat¬ 
tice models with pair interaction between the first neigh¬ 
boring molecules, with interactions spanning two fl3| and 
three lattice sites El- While in Ref. [l3] both van der 
Waals and hydrogen bond like interactions were used, 
resulting in a line of temperature of maximum density 
associated to a ground state phase transition (GSPT), in 
Ref. El ^ was proposed a core-softened fluid with pair 
interactions up to three sites, resulting in two tempera¬ 
ture of maximum density lines associated to two GSPT. 
Besides obtaining exact results, in the latter work we 
used an analytical approximation in the neighborhood of 
the critical point to obtain a simple expression for Gibbs 
free energy, and used it to mathematically study the re¬ 
lation between anomalous density behavior and GSPT. 


In this work we proceed on this direction by investi¬ 
gating the repulsive ID lattice gas, which is even simpler 
than our previous models and presents waterlike anoma¬ 
lies in density, thermodynamic response functions and 
self diffusion constant. The model was studied through 
transfer matrix technique, the Takahashi method (within 
a two state approximation, as will be discussed latter), 
and Monte Carlo simulations. With the results obtained 
from these techniques a connection between temperature 
of maximum density and GSTP was found as in a pre¬ 
vious work with more complex models El- In addition, 
it was also found that GSPT does present a residual en¬ 
tropy, due to phase mixing, and it is shown that this 
property is fundamental in determining waterlike anoma- 
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lies for the model considered here. Finally, a comparison 
between regions with density and diffusion anomaly in¬ 
dicated that this model presents so called hierarchy of 
anomalies [Tst j . 


The eigenvalue A is related to the fugacity 2 = e^ via 
characteristic function of V, 

(az — A) (1 — A) — z = 0, (5) 


II. THERMODYNAMICS 
A. Model and Ground State 

By defining the presence of a particle on site k with 
an occupation variable rjk = 1 (and absence with r)k = 0) 
the effective Hamiltonian of the ID repulsive lattice gas 
in the grand canonical ensemble becomes: 

n = ^i+i - Wi, (i) 

i i 

where e > 0 is the strength of the (repulsive) interaction 
and /i is the chemical potential. Note that this Hamilto¬ 
nian can be seen as a special case from our previous model 
for liquid water 0 by setting e hb = 0 and e v d w = —£■ 

The ground state of the unidimensional repulsive lat¬ 
tice gas is composed by three phases: a gas in which the 
lattice is empty, a dense fluid (DF), presenting a com¬ 
pletely filled lattice, and a half-filled lattice with each 
particle separated by holes, here denominated as soft¬ 
ened fluid (SF)[h|- In a lattice with L sites separated 
by a distance Z, periodic boundary conditions, and N 
particles, the enthalpy per particle of these phases are 
Zigas = 0, /idf = e + PI and h SF = 2 PI, with P meaning 
pressure. Considering this, the SF is stable in the ground 
state for pressures 0 < P < P c with P c l = e. 

Two ground state phase transitions are present at null 
temperature and were calculated by equating the en¬ 
thalpy per particle at different phases: a G-SF at P = 0 
and a SF-DF transition at P c , which will be called “sec¬ 
ond critical point” as in previous works 0E1II3. 


B. Transfer Matrix Technique 

Thermodynamics of the model proposed here is ob¬ 
tained using transfer matrix technique and will be de¬ 
scribed shortly. Other derivations using this technique 
can be found elsewhere [Is|- We start from the grand 
canonical partition function as a trace of a matrix 

S(T,L, M ) = ^e-^ = Tr{P L }, (2) 

v 

where /3 = 1/kpT, fj = { 771 ,..., t]n}, and the elements of 
V are given by 

TV =e- /3 ( ew '" H. (3) 

In the thermodynamic limit only the largest eigenvalue 
of V will contribute to the partition function: 

H(T,L,/i) = \ L =e f)plL . (4) 


with a = e -/3e . 

Next we introduce reduced variables t = kpT/e, p = 
Pl/e, g = /r/e, v = V/Nl and s = S/Nkp • Since the 
equation of state (0 connects the Gibbs free energy g to 
thermodynamic variables t and p , it is possible to obtain 
the ‘volume’ v and entropy s per particle as 



With these definitions, density can be calculated as 


1 (1 — A) [a (1 — A) — 1] 

v a(l — A ) 2 — 1 + 2A 


( 6 ) 


and entropy becomes 


s 


, r \ , P , za i0- ~ X )+P X \ ~P X2 
lll(j) + 7 + rt [a(l — A) — 1| 


( 7 ) 


Thermodynamic response functions, such as thermal ex¬ 
pansion coefficient a, isothermal compressibility Zct, and 
isobaric heat capacity cp are calculated from standard 
definitions [31 which, within the convention adopted 
here, reads 



( 8 a) 

( 8 b) 

( 8 c) 


Explicit (and exact) expressions for these functions are 
too lengthy to be reproduced here. The temperature of 
maximum density line was calculated from () 8 a|l by nu¬ 
merically solving av = 0 . 


C. Two states approximation 

To obtain further theoretical insight about this system, 
it is interesting to investigate its entropy in the neigh¬ 
borhood of the second critical point using a two state 
approximation based on a lattice version of the Taka- 
hashi method, as employed in Ref. 0 . In that work the 
Gibbs free energy per particle was approximated near 
the GSPT by realizing that on this region the proba¬ 
bility of finding local states different from those of the 
coexisting structures was vanishing small. In the cur¬ 
rent model these states would correspond to softened and 
dense fluids, whose ground state local molecular spacing 
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FIG. 1. Location of temperature of maximum density line 
( Tmd ) in the p vs. t phase diagram. The ground state phase 
transition between the softened fluid and dense fluid is indi¬ 
cated with a triangle, while its pressure value is marked with a 
dotted line. An example of the mixed SF/DF state is shown 
along the dotted line, used to indicate the critical pressure 

pm Pc = l. 


are rs = 2 m = 21. Thus, the Gibbs free energy (non- 
reduced units) becomes approximately: 


G/N « h c ~k B T In 


e ~P (P-PcVs + e -P(P~Pc)r D 


(9) 


where h c is a constant identical to the local microscopic 
enthalpy at the GSPT (see Ref. ,3l)- Within the no¬ 
tation adopted here, the Gibbs free energy on Eq. © 
reads 


h c 3 

q = -1-An — t In 

y e 2 


2 corf, (^) 


( 10 ) 
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FIG. 2. The temperature dependence of (a) density and 
(b) entropy for pressures equal, slightly below and above the 
ground state phase transition at p = p c = 1. 


with A p = p — p c . From the latter, molecular volume 
and entropy becomes 


3 1 

v = -tanh 

2 2 



( 11 ) 


and 


s = In 


2 cosh 




( 12 ) 


It is possible to see that by approaching the GSPT at 
constant pressure, i.e., (t —>• 0,p = p c = 1), one obtains 
the residual entropy s —» In 2, discussed before. For any 
other route used to approach the ground state, entropy 
results in a null value, indicating that in this model resid¬ 
ual entropy only occurs exactly at the GSPT. Molecular 
volume also depends on the path along which criticality 
is approached, resulting in 


lim lim v = 

P^rPc t—>o 


2 , 

1. 


lim lim v = —. 

t—yOp—>p c 2 


P -t Pc 

P^-Pt 


(13a) 

(13b) 


It is important to investigate the critical behavior of the 
system near the GSPT, by calculating analytical expres¬ 
sions for thermodynamic response functions (j8al) - (|8c|> . 
within the two-states approximation adopted here. It is 
convenient to define f(x ) = 1 — tanh 2 x, with x = Ap/2f, 
such that 


av = A p 1 x 2 /(x), (14a) 

ktv=^-, (14b) 

Cp = x 2 f{x). (14c) 

Since 0 < /(x) < 1, thermodynamic response functions 
diverge toward criticality as kt ~ £ -1 , a ~ A p~ 1 x 2 and 
c p ~ x 2 . The low temperature behavior on response func¬ 
tions is consistent with waterlike behavior and the critical 
exponents for k t and c p are identical to those found in 
a continuous one dimensional core softened model pro¬ 
posed by Sadr-Lahijany et al E2 in the same context. 
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FIG. 3. Entropy as a function of pressure at a fixed tem¬ 
peratures from exact (continuous and dashed lines) and two 
state approximation (circles). Note that entropy presents a 
maximum value at low temperatures. 



FIG. 4. Thermal expansion coefficient, a, as a function of 
pressure at fixed temperatures. At low temperatures, an oscil¬ 
latory behavior is observed in the neighborhood of the ground 
state critical pressure p c = 1. 


D. Fluid Structure 

We start discussing our model’s liquid structure from 
its exact pressure vs. temperature phase diagram, shown 
in Fig. [1] The model presents a GSPT between softened 
and dense phases at p = p c = 1, from which a tempera¬ 
ture of maximum density line, Tmd , emerges separating 
a region with density anomaly (a < 0) from a region 
with normal density behavior (a > 0). This line was ob¬ 
tained by numerically solving the equation a = 0, and 
was found to extend up to a temperature t = 0.279(1), 
and proceed retracing to lower temperatures ending ex¬ 
actly at pressure pi ow = 1/2. An explanation for this 


precise number can be obtained is as follows: at null tem¬ 
perature, the entropy gain for contracting or expanding 
the softened fluid by a single site will be same (a lattice 
defect), thus the system contracts or expands through a 
competition between fluctuations in energy and volume, 
which mathematically corresponds to free energy fluctu¬ 
ations of 8g con = e — PI (contraction) or Sg exp = +PI 
(expansion). By equating both variations it follows that 
Plow — 1/2- 

Another interesting feature of this phase diagram is the 
mixed SF/DF state, which occurs near the GSPT. This 
can be observed in Fig. [2] (a), where density is shown as 
a function of temperature for the critical pressure p c = 1 
and values slightly below and above p c . Density anomaly 
is evident for p < p c while density behaves normally for 
p > p c . Exactly at p = p c = 1 density is kept constant 
up to t fa 0.15. Density behaves this way because a oc 
(bp/dT)p 0 near the Tmd line. 

The entropy increase with temperature is depicted in 
Fig. [2] (b), for the pressures investigated in Fig. [2] (a). 
From this data it is evident that a residual entropy occurs 
exactly at p = p c = 1, while for other pressures ground 
state entropy start to increase from a null value. It is 
well known that, due to the Maxwell relation 

-(S)„-©„ » 

a density maximum at fixed pressure is mathematically 
equivalent to an entropy maximum as function of pres¬ 
sure at fixed temperature. Thus, in Fig. [3j entropy is 
investigated as a function of pressure (at fixed tempera¬ 
tures) using exact and approximate expressions obtained 
previously, Eqs. 0 and (1121) . It is important to note 
accuracy of the two states approximation near the tran¬ 
sition. The relevance of a residual entropy for the appear¬ 
ance of a Tmd line is evident since a region of anomalous 
entropy increase appears at low temperature in the neigh¬ 
borhood of the SF/DF transition. What is interesting in 
this behavior is that entropy anomaly appears naturally 
due to residual entropy, whose origin is the two phase 
mixture in the GSPT. In this way, anomalous entropy 
behavior (and anomalous density) is intrinsically associ¬ 
ated to phase transition in our ID model since mixing 
and critical behavior cannot be dissociated when contin¬ 
uous and discontinuous transitions are collapsed in the 
same point. 

Residual entropy also explains the oscillatory behavior 
observed on thermal expansion coefficient (a) as a func¬ 
tion of pressure, near the critical transition, Fig. [4] This 
behavior was observed previously 00 and its rela¬ 
tion with entropy anomaly is simple: below the critical 
pressure a must be negative due to the fast increase in 
entropy towards its residual value In 2, and above p c en¬ 
tropy must fast decrease to a low value, since it is null at 
t = 0. 

It’s worth mentioning that two states approximation 
also predicts that isothermal compressibility and con¬ 
stant pressure heat capacity increase while lowering tem- 
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FIG. 5. Self-diffusion constant as a function of density at 
fixed temperature. 


FIG. 6. The location of the maximum (circles) and minimum 
(square) values of the self-diffusion constant are compared to 
the line of temperature of maximum density (line) in the t vs. 
p phase diagram. 


perature, at pressures near the critical one, through 
Eqs. (114bl) and (114cf> . Both response functions must 
present maximum values as a function of temperature, 
for pressures near the critical one. These results are not 
shown here but they can be rationalized as follows: the 
increase in isothermal compressibility should be expected 
near the GSPT because while increasing pressure in that 
region the molecular volume is abruptly changing from 
typically softened fluid to typically dense fluid. On the 
other hand, the existence of a maximum as a function 
of temperature in c p , at pressures near the critical value, 
can be predicted from Eq. (|8c]> and Fig. [2] (b). Its ex¬ 
planation proceeds as follows: let us consider a system 
in the softened fluid, at t = 0 and p < p c , where all pair 
of particles are separated by a hole. In this state, while 
increasing the temperature entropy will preferentially be 
increased by allowing neighboring particles to contract. 
As temperature increases further, an even higher number 
of particles pairs are allowed to fluctuate between soft¬ 
ened and dense fluid states, but the nearest the system 
is to the critical pressure, the fastest is the approach to 
an entropic value of In 2, since the enthalpy difference 
between the two states is proportional to A p. At higher 
temperatures this mixing mechanism will compete with 
translational diffusion, i.e., volume expansion, there must 
be a maximum in ( ds/dt) p between those temperatures 
were each mechanism is dominating the entropic behav¬ 
ior of the system. 


III. MONTE CARLO 

Standard NVT Monte Carlo simulations with 
Metropolis algorithm were used to calculate the self 
diffusion constant of our model as a function of density 
and temperature . The system dynamics proceeds 

as follows: particles seated on site k were allowed to 


jump to nearest neighbor (k ± 1) and next nearest 
neighbor sites (fc±2). Jumping to next nearest neighbor 
sites is necessary for particles to diffuse in a ID lattice, 
because there is no way for them to turn around their 
neighboring particles (the same strategy was used in our 
previous paper [l3|. System size was equal to L = 10 3 
and simulation times were equal to or higher than 10' 
Monte Carlo steps after an equilibrium time of 10 6 
Monte Carlo steps. 

As in Refs. [13|, |22|, |23|], self diffusion constant was ob¬ 
tained by linearly adjusting the mean square displace¬ 
ment of particles to the Einstein-Smoluchowski equation 
for the Brownian motion in one dimension, 

1 W \ 

-J2iX k (r)~X k ( 0)] 2 \ =2 D(t,p)T (16) 

k=l / u 

where X k (r) is the position of particle k at time r and 
1Z > 5 is the number of different initial conditions used 
to average square displacements. 

Fig. [5] exhibits D{p 1 t) for different isotherms. At 
higher temperatures, as illustrated for t = 0.7, the dif¬ 
fusion constant is ‘normal’ and monotonically decreases 
with density. For t < 0.55 minimum and maximum 
diffusion appears, delimiting a region where diffusion 
anomalously increase with density, as observed experi¬ 
mentally 0. 

In liquid water, density anomaly is associated to dif¬ 
fusion anomaly through a hierarchy of anomalies, with 
the former enclosing the latter in the phase diagram 0- 
To investigate this hierarchy diffusion isotherms were ad¬ 
justed to a fourth degree polynomial in the neighborhood 
of extremum temperatures. At low temperatures, the 
location of maximum and minimum diffusion were esti¬ 
mated using six points around each simulated extremum, 
and the best estimate for their values were taken from the 
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analytically adjusted function. On intermediate temper¬ 
atures, where maximum and minimum where too close, 
a single function and a higher number of points were 
used to adjust diffusion isotherms. The result is shown 
in Fig. 1 where is possible to see that an strict hierarchy 
of anomalies is found, with the same order observed in 
atomically detailed models of liquid water [l5j . 

IV. FINAL DISCUSSIONS AND CONCLUSIONS 

We investigated the exact thermodynamics of the 
one dimensional repulsive lattice gas model using trans¬ 
fer matrix technique. In the ground state, the model 
presents a phase transition between a softened fluid and 
a dense fluid, which is characterized by a mixture be¬ 
tween two states. The main consequence of this mixed 
state is a residual entropy in a single point in the pressure 
vs. temperature phase diagram. Using thermodynamic 
relations, it was argued that the presence of a residual en¬ 
tropy in a single point is consistent with a temperature 
of maximum density line (a = 0) emanating from the 
ground state phase transition. Considering this it is pos¬ 
sible to state that ground state phase transition, residual 
entropy and density anomalies are thermodynamic prop¬ 
erties of our model which are intrinsically related and 
cannot be dissociated from each other. 

Further discussion about ground state phase transi¬ 
tions in one dimensional models is relevant at this point. 
First of all, the two transitions observed in our model can 
be seen as remanent from phase transitions on higher di¬ 
mensional systems collapsing to T —> 0 as dimensionality 
is reduced to d —> 1. This collapse is also evident on 
the analytical solution of a lattice gas inside the Bethe 
lattice: in this case, while continuously reducing the co¬ 
ordination number the liquid-gas phase transition con¬ 
tinuously shrinks to T —> 0 [l]a, (25( • This is a relevant 
issue for one dimensional systems because the observed 
ground state phase transitions incorporate elements from 
both continuous and discontinuous phase transitions [26| . 
While usually a discontinuous transition allows for any 
mixtures between two coexisting phases, when a contin¬ 
uous transition is collapsed on it, any mixture except 
for that maximizing the entropy is forbidden. Thus, 
for fluid-fluid transition entropy is maximized for equal 
distribution of randomly mixed interactions compatible 
with softened and dense fluids, resulting in an entropy 
per particle s = S/N = kb In2 (as in a coin toss game). 

The viability of using a simple two states approxima¬ 
tion to accurately describe the system in the neighbor¬ 


hood of the ground state phase transition is another issue 
deserving attention, since there is a long tradition on us¬ 
ing two state models to describe liquid water [H HgJ ■ 
Recently, a modified regular solution model within a two 
solvent approach was adjusted to computer simulation 
data from ST2 model, being accurately used to describe 
computational data above the second critical point, in a 
region called Widom line [29J. This raises the possibil¬ 
ity of connecting the current approach to more complex, 
higher dimensional and off-lattice models. Nevertheless, 
it should be emphasized that the two states character of 
the approximation employed here is conceptually distinct 
from the ones implemented in Ref. [2f| , where two states 
refers to two fluids with distinct free energies. Our ap¬ 
proach resembles more the approximation implemented 
in Ref. [§}, in the sense that two states are localized 
microscopic states with well defined local energies and 
entropies, which are useful for truncating the partition 
function in some kind of expansion. Certainly it should 
be interesting to link these different views, but this is 
outside the scope of the current work. 

To finish, it should be mentioned that in a previous 
work on lattice models with waterlike behavior [TiH we 
found an oscillatory behavior on thermal expansion coef¬ 
ficient a associated to density anomaly and ground state 
phase transition. More recently this effect was shown to 
be general, being a feature of any ground state phase 
transition between different fluid structures in one di¬ 
mensional lattice models with interactions restricted to 
first neighboring particles and, inspired by these models, 
it was possible to design a three dimensional spherically 
symmetric pair potentials with two liquid-liquid phase 
transitions, both of them associated to temperature of 
maximum density lines 0 - The oscillatory behavior on 
a , as well as the anomalous behavior on isothermal com¬ 
pressibility and constant pressure heat capacity, can be 
explained in terms of a residual entropy using thermody¬ 
namical and statistical arguments, while for the current 
model the the relation between GSPT, residual entropies 
and waterlike behavior can be mathematically described 
within the two states approximation. 

Obviously we are tempted to ask whether such entropic 
effects couldn’t play an important rule in more complex 
systems, such as higher dimensional and off-lattice mod¬ 
els. Current work is being done on this direction through 
molecular dynamics simulations of the 3D core softened 
models, as in Ref. fl4j . 
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